import sys
sys.path.append('c:\\Users\\karen\\Dropbox\\DataQuality\\libsvm-3.1\\python')
sys.path.append('c:\\Users\\karen\\Dropbox\\DataQuality\\libsvm-3.1\\windows')

import svm
from svm import *
from multiprocessing import Pool
from numpy import *
import matplotlib
matplotlib.use("Agg")
from pylab import *

def grid_regress (nurange,crange,gammarange,output,input,nf,timeout,p):
	params = []
	for nu in nurange:
		for c in crange:
			for gamma in gammarange:
				params.append([nu,c,gamma,nf,output,input])

	N = len(params)
	results = p.imap(do_one_cv,params)

	del params
	maxscc = 0
	minerror = inf
	ind = 0
	optgamma = 0
	optc = 0
	optnu = 0
	mses = {}
	
	for nu in nurange:
		for c in crange:
			for gamma in gammarange:
				if ind % (N * 0.05) == 0:
					print '====== %%%g completed ======' % (float(ind)/N * 100)
				ind += 1
				try:
					mse,scc = results.next(timeout=timeout)
				except:
#					print '(%g,%g,%g)' % (nu,c,gamma),
#					print '* ',
					mse = inf
					scc = 0
				mses[nu,c,gamma] = mse
				if scc > maxscc or scc == maxscc and gamma == optgamma and c < optc and nu < optnu:
					print scc, mse,gamma, c, nu
					minerror = mse
					maxscc = scc
					optc = c
					optgamma = gamma
					optnu = nu
			


	del results
	return [maxscc,minerror,optnu,optc,optgamma,mses]


def grid_classify (crange,gammarange,output,input,nf,useprob,timeout,p):
	params = []
	for c in crange:
		for gamma in gammarange:
			params.append([c,gamma,nf,output,input,useprob])

	N = len(params)
	results = p.imap(do_one_cv_classify,params)

	
	del params
	maxauc = 0
	ind = 0
	optgamma = 0
	optc = 0
	optbias = 0.5
	maxoptacc = 0
	optacc = 0
	bias = 0
#	accs = {}
		
	for c in crange:
		for gamma in gammarange:
			if ind % (N * 0.05) == 0:
				print '====== %%%g completed ======' % (float(ind)/N * 100)
			ind += 1
			try:
				auc,optacc,bias = results.next(timeout=timeout)
			except:
#					print '(%g,%g,%g)' % (nu,c,gamma),
#					print '* ',
				auc = 0
#			accs[c,gamma] = acc
#			print auc,optacc,bias
			if auc > maxauc or auc == maxauc and gamma == optgamma and c < optc:
				print 'Max AUC = %g (gamma = %g, c = %g)' % (auc,gamma, c)
				maxauc = auc
				maxoptacc = optacc
				optbias = bias
				optc = c
				optgamma = gamma


	del results
	return [maxauc,maxoptacc,optbias,optc,optgamma]#,accs]
	
	
	
	
def do_one_cv(theinput):
	nu = theinput[0]
	c = theinput[1]
	gamma = theinput[2]
	nf = theinput[3]
	output = theinput[4]
	input = theinput[5]

	param = svm.svm_parameter('-s %d -t %d -n %g -c %g -g %g' % (svm.NU_SVR,svm.RBF,nu,c,gamma))

	prob = svm.svm_problem(output, input)
	target = (c_double * prob.l)()
	
	libsvm.svm_cross_validation(prob, param, nf, target)	
	MSE,SCC = evaluations(prob.y[:prob.l],target[:prob.l])
	del target
	return MSE,SCC

	
def do_one_cv_classify(theinput):
	c = theinput[0]
	gamma = theinput[1]
	nf = theinput[2]
	output = theinput[3]
	input = theinput[4]
	useprob = theinput[5]

	param = svm.svm_parameter('-c %g -g %g -b %d' % (c,gamma,int(useprob)))

	prob = svm.svm_problem(output, input)
	target = (c_double * prob.l)()
	posclass = output[0]
	libsvm.svm_cross_validation(prob, param, nf, target)
	ys = prob.y[:prob.l]
	db = array([[ys[i],target[i]] for i in range(prob.l)])
	
	del target
	
	neg = len([x for x in ys if x != posclass])
	pos = prob.l-neg;
	auc,topacc,optbias,top_tps_bias,top_fps = calc_AUC(db,neg,pos,posclass)#,min_fps_in_tp=[],do_plot=False,pagetitle='',pdfpage=0)
	
	return auc,topacc,optbias


	
def evaluations(ty, pv):

	if len(ty) != len(pv):
		raise ValueError("len(ty) must equal to len(pv)")
	sumv = sumy = sumvv = sumyy = sumvy = 0
	total_error = 0

	for v, y in zip(pv, ty):
		total_error += (v-y)*(v-y)
		sumv += v
		sumy += y
		sumvv += v*v
		sumyy += y*y
		sumvy += v*y 
	l = len(ty)
	MSE = total_error/l
	try:
		SCC = ((l*sumvy-sumv*sumy)*(l*sumvy-sumv*sumy))/((l*sumvv-sumv*sumv)*(l*sumyy-sumy*sumy))
	except:
		SCC = float('nan')
	return (MSE, SCC)


def evaluations_classify(ty, pv,posclass=-1,bias=0.5):

	if len(ty) != len(pv):
		raise ValueError("len(ty) must equal to len(pv)")
	sumv = sumy = sumvv = sumyy = sumvy = 0
	total_correct = 0
	
	confusionmatrix = zeros([2,2]);
	
	for v, y in zip(pv, ty):
		if y == posclass and v>bias: 
			confusionmatrix[0,0] += 1
			total_correct += 1
		elif y==posclass and v<bias:
			confusionmatrix[0,1] += 1
		elif y!=posclass and v>bias:
			confusionmatrix[1,0] += 1
		elif y!=posclass and v<bias:
			confusionmatrix[1,1] += 1
			total_correct += 1
			
	l = len(ty)
	ACC = 100.0*total_correct/l
	return ACC,confusionmatrix
	

def calc_AUC(db,neg,pos,posval,min_fps_in_tp=[],do_plot=False,pdfpage=0,pagetitle = ''):

	allNs = range(size(db,0))
	shuffle(allNs)

	db = db[allNs,:]
#	print db
	db1 = copy(db)
	db = db[argsort(db[:,1]),:]
	
#	print db
	
	xy_arr = [[1,1]]
	topacc = -1
	tp,fp = pos, neg
	#assure float division
	N = neg+pos
#	print N
	accs = [[(tp+(neg-fp))*1.0/N,0]]

	defailt_bias_ind = -1

	prev_tp = 1.0
	top_fps = []
	top_tps_bias = []
	min_fps_ind = -1

	for i in range(size(db,0)):
		if db[i,0]==posval:                 #positive
			tp-=1.0
		else:
			fp-=1.0


		for min_fp_in_tp in min_fps_in_tp:
			if tp/pos < min_fp_in_tp and prev_tp >= min_fp_in_tp:
				if i == 0:
					top_tps_bias.append(0)
					top_fps.append(fp/neg)
					min_fps_ind = size(db,0)-i
				else:
					top_tps_bias.append(db[i-1,1])
					top_fps.append(fp/neg)
					min_fps_ind = size(db,0)-i

#               print tp/pos,fp/neg,prev_tp,min_fp_in_tp,min_fps_ind,i


		prev_tp = tp/pos

		xy_arr.insert(0,[fp/neg,tp/pos])
		accs.insert(0,[(tp+(neg-fp))*1.0/N,db[i,1]])

		if db[i,1] > 0.5 and defailt_bias_ind == -1:
			defailt_bias_ind = size(db,0)-i

	topacc = max([x[0] for x in accs])
	topacc_bias = [i for i in range(len(accs)) if accs[i][0] == topacc]


	auc = 0.
	prev_x = 0
	for x,y in xy_arr:
		if x != prev_x:
			auc += (x - prev_x) * y
			prev_x = x


#	print xy_arr
	if do_plot:
		plot([x[0] for x in xy_arr],[x[1] for x in xy_arr])
		plot(xy_arr[topacc_bias[0]][0],xy_arr[topacc_bias[0]][1],'or')
		plot(xy_arr[defailt_bias_ind][0],xy_arr[defailt_bias_ind][1],'ob')
		plot(xy_arr[min_fps_ind][0],xy_arr[min_fps_ind][1],'og')


		xlabel("False Positive Rate")
		ylabel("True Positive Rate")
		title("%s (AUC = %.4f; ACC = %.4f)" % (pagetitle,auc,topacc))
		grid()
		pdfpage.savefig()
#		'%s_ROC.pdf' % (outputfn),format='pdf')
		close()

#	print auc,topacc,accs
	
	return [auc,topacc,accs[topacc_bias[0]][1],top_tps_bias,top_fps]